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On the Pulse Intensity Modulation of PSR B0823+26 
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ABSTRACT 

We investigate the radio emission behaviour of PSR B0823+26, a pulsar which is known to 
undergo pulse nulling, using an 153-d intensive sequence of observations. The pulsar is found 
to exhibit both short (~ min) and unusually long-term (~ hours or more) nulls, which not 
only suggest that the source possesses a distribution of nulling timescales, but that it may 
also provide a link between conventional nulling pulsars and longer-term intermittent pul- 
sars. Despite seeing evidence for periodicities in the pulsar radio emission, we are uncertain 
whether they are intrinsic to the source, due to the influence of observation sampling on the 
periodicity analysis performed. Remarkably, we find evidence to suggest that the pulsar may 
undergo pre-ignition periods of 'emission flickering', that is rapid changes between radio-on 
(active) and -off (null) emission states, before transitioning to a steady radio-emitting phase. 
We find no direct evidence to indicate that the object exhibits any change in spin-down rate 
between its radio-on and -off emission modes. We do, however, place an upper limit on this 
variation to be < 6 % from simulations. This indicates that emission cessation in pulsars does 
not necessarily lead to large changes in spin-down rate. Moreover, we show that such changes 
in spin-down rate will not be discernible in the majority of objects which exhibit short-term 
(< 1 d) emission cessation. In light of this, we predict that many pulsars could exhibit simi- 
lar magnetospheric and emission properties to PSR B0823+26, but which have not yet been 
observed. 

Key words: methods: data analysis - pulsars: individual: PSR B0823+26 - pulsars: general. 



1 INTRODUCTION 

PSR B0823+26 was one of the first pulsars to be discovered 
dCraftetal.li968h and, as such, has been studied for over 40 yr. De- 
spite its typical spin parameters (i.e. rotational period P ~ 0.53 s, 
period-derivative P ~ 1.71 x 10 -15 and magnetic field strength 
B ~ 0.96 TG), this pulsar is by no means ordinary. Among its 
most salient features, are its int er-pulse and post-cursor emission 
components l lBacker et alj 1973b . which are typically obse rved at a 
few percent of its main-pulse peak intensity dRankml 19861 : see also 
Fig.[2j. These features are quite rare among the pulsar population, 
particularly the inter-p ulse emission which is only observed among 
a small subset (~ 3 %: lMaciesiak et ai]|201 lh of the normal pulsar 
population. P SR B0823+26 is al so found to exhibit pulsed, soft X- 
ray emission ( iBecker et ai1l2004l) . which makes it only one out of 
nine old pulsars (1 Myr < r < 20 Myr) that produce high-energy 



emission 



(e.g. Becker et"al] |2004| ; iTepedelenlioglu & Ogelmanl 
l2005tlLietalj2008l) . 

In addition to these rare emission features, the object ex- 
hibits abrupt cessation and re-activation of its radio emission 



I Hesse & Wielebinskill9"74l : lRitchingsll976l : lRathnasree & Rankin 



1 1995b . that is commonly referred to as 'pulse nulling' backed 
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Il970h . This phenomenon, which is also observed in numerous 
other pulsars, a ffects all components of emission dBackerl 1 1970l 
iRitchingsl fl976h and can be thought of as an extreme mani- 
festation of mode changing, that is a variation between active 
(radio-on, hereafter) and quiescent (radio-off or null, hereafter) 
emission modes. Through studying the emission behaviour of 
PSR B0823+ 26 over three observing runs (< 2 h of observa- 
tions in total). iRathnasree & Ran kin (1995) infer a nulling frac- 
tion (NF) — the fraction of pulses which are radio quiet — that 
is 6.4 ± 0.8%. 



While th e documented n ulls in PSR B0823+26 last a few 
pulse periods dRitchingslll976h . pulse nulling can be observed over 
a much wider range of timescales in pulsars; that is, from just 
one or two pulses to many days (e. g. lRankinlll986l : iBiggslI 19921: 
iKramer et alj2006l : IWang et alj2007h . Accordingly, the NFs of pul- 
sars can also range from less than 1 % up to, in excess of, 95 % 
dPeich et al]|l986t IWang et alj|2007l) . Due to this diversity in the 
observed nulling statistics, and incomplete sample of known ob- 
jects which exhibit temporary radio emission failure — < 200 
nulling (transient) pulsars are currently known, and even less have 
been studied in detail — not much is known about the 'typical' 
properties of nulling pulsars (e.g. P, P or magnetic inclination an- 
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gle afl nor why exactly they undergo such extreme changes in 
their radio emission mechanism. 

It has been suggested that 'magnetospheric-state switching' 
could be the underlying mec hanism responsible for mode-changing 
and nulling in p ulsars (e.g. iBartel et al.|[l982l ; IContopoulosll2005l ; 
iTimokhinlHoTa) . This process refers to alterations to the global 
current distribution in a neutron star magnetosphere which, in 
turn, are thought to cause moding between different emission 
states. These global, magnetospheric reconfigurations are also pro- 
posed to explain the (quasi-)periodic timing signatures in pul- 
sar residuals through correlated changes in spin-d own rate {u), 
such as those in the intermittent pulsar B 1931+24 ^Kramer et al.1 
l200dlHobbs et al.ll2010l : iLvne et al.ll2010h . However, not much is 
known about how such magnetospheric alterations are triggered, 
nor what their physical timescales or periodicities should be; there 
are a number of different models that a ttempt to explain th is phe- 
nomenon, e.g. non-radial oscillatio ns dRosen et al.l 1201 it) , aster- 
oid belts dCordes & Shannonl [2008h . precessional torques Jjonej 
|2012|). surface tem perature variations in the pol ar gap region 
dZhang et alJI 1997 ) and magnetic field instabil i ties {Geppert et al .1 



l2003l ; IUrpin & Gill2004lRheinhardt et alj2004lWang et al.l2007l) . 
but none are conclusive. 

With the above in mind, we present the analysis of recent high- 
cadence observations of PSR B0823+26, which were stimulated 
by the discovery of longer than previously recorded nulls during 
regular timing observations with the Lovell Telescope at Jodrell 
Bank. These data represent the most comprehensive study of the 
intermittent behaviour of this pulsar for the first time. This sub- 
sequently enables the comparison between long-term intermittent 
pulsars and more 'conventional' nulling pulsars, as well as insight 
into the mechanism(s) which govern radio emission modulation in 
such objects. In Section [2] we describe the observations, followed 
by an overview of the emission variability of the source in Sec- 
tion[3] In Section[4] we review the timing behaviour of the neutron 
star, and the simulation tool used to model its rotational properties. 
Finally, we discuss the implications of the results in Section [5] and 
present our conclusions in Section[6] 



2 OBSERVATIONS 

The observations of PSR B0823+26, presented here, were obtained 
with the Lovell Telescope over a period of approximately 153 d 
(1 January 2009 to 3 June 2009). These data were obtained using 
the Analogue Filter Bank (AFB) and Digital Filter Bank (DFBfj 
back-ends using a 1400 MHz receiver, with an average cadence of 
approximately 4 observations per day (see TableQ]for details). 

In addition to this typical daily monitoring, we carried out 
three observing runs with the AFB which spanned the entire time 
the source is above the horizon at Jodrell Bank (see Tablef2]for de- 
tails). These data were obtained to compare the long-term variation 
with possible short timescale modulation. As a result, there are sig- 
nificantly more observations with this back-end compared with the 



1 While the NF of a pulsar i s found to be we akly correlated with 
characteristic age iRankirJ [l98el IWang et al.1 12007|) . no significant cor- 
relation has been found with other basic pulsar properties or geometry 
iRankin & Wrighl2007l) . 

2 The DFB was commissioned approximately 20 d after the start of the 
AFB observations. 



Table 1. System characteristics of the observations of PSR B0823+26 
described here. The bandwidth capability of the DFB was increased on 
1 April 2009, i.e. MJD ~ 54922.1. Therefore, values quoted for the DFB 
back-end are with respect to times prior to (DFB prc ) and post (DFBp OS t) 
MJD ~ 54922.1 respectively. 



System Property 


AFB 


DFBp le 


DFBp OS t 


Time span of observations (d) 


153.3 


77.3 


60.9 


Total number of observations 


1274 


284 


281 


Typical observation duration (min) 


6 


6 


6 


Average observation cadence (d -1 ) 


8.3 


3.7 


4.6 


Typical sky frequency (MHz) 


1402 


1382 


1374 


Typical observing bandwidth (MHz) 


32 


113 


128 


Typical channel bandwidth (MHz) 


1 


0.25 


0.25 



Table 2. The observation properties for the three continuous observing in- 
tervals during January 2009. Each continuous observing run has a length 
T oas , for which there are A r or)S separate observations. 



Start Epoch (MJD) 


T obs (d) 


Nobs 


54837.7215 


0.6390 


152 


54854.6844 


0.6623 


155 


54857.7146 


0.6067 


145 



DFB. Despite the deficit in the number of observations, the data ob- 
tained with the DFB is complementary to the emission modulation 
study due to the better sensitivity and wider bandwidth. 



3 EMISSION VARIABILITY 
3.1 Nulling activity and flux limits 

The denser coverage provided by the AFB observations meant that 
they were used for characterising the overall emission modulation 
properties of PSR B0823+26. To determine whether the pulsar was 
radio-on or -off, we used the average profiles formed over the entire 
bandwidth for each approximately 6-min long observation. This 
was done by visual inspection and resulted in a time-series of one- 
bit data corresponding to the 'radio activity' of the pulsar i.e. l's for 
observations with detectable emission and 0's for observations with 
non-detectable emission (i.e. radio-off states fj To complement the 
visual inspection, we also made use of the timing model to con- 
firm detections. A time-of-arrival (TOA) was calculated for each 
observation and was compared to the known timing model. Those 
in good agreement were confirmed as detections. 

FigureQ]shows an example of an observation where the pulsar 
is observed to transition between the radio-on and radio-off phases. 
The transition timescales between the emission phases for this pul- 
sar are similar to those seen in PSR B 193 1+24; that is, of the or- 
der of seconds or less. Such sharp discontinuities in pu lse intensity 
are in consistent with interstellar scintillation (see, e.g. JWang et al.l 
1 20051 and references therein). Therefore, these transitions between 
emission phases are considered to be intrinsic to the pulsar. 



3 We note that only the nulling activity of the main-pulse component is 
considered here, as the 6-min observations do not provide enough sensitiv- 
ity to characterise the inter-pulse or post-cursor emission properties. 



On the Pulse Intensity Modulation ofPSR B0823+26 3 




50 100 150 200 250 300 3 50 



Pulse Longitude (deg) 

Figure 1. Consecutive pulse profile sub-integrations for PSR B0823+26 
(from bottom to top), which were obtained during one observation on 
10 May 2009. The pulse intensity of each sub-integration is normalised 
to one and is offset from the next profile for clarity. The pulsar is de- 
tectable in the first four sub-integrations (~ 60 s each), after which it 
abruptly 'switches off'. This transition occurs sometime during the fourth 
sub-integration over a timescale of less than one minute. 

Due to the greater sensitivity, we use the DFB data to pro- 
vide limits on the average pulse flux density during the sep- 
arate phases of emission. By considering an integrated profile, 
formed from an observation of length T, with an equivalent pulse 
width W cq and signal-to-noise ratio SNR, the mean flux den- 
sity can be estimated vi a the modified radiometer equation (e.g. 
lLorimer & Kramej|2005h : 

/3SNRr sys / W cq 

Here, /3 ~ 1 is the digitisation factor, G ~ 1 JyK -1 is the tele- 
scope gain, T sys ~ 35 K is the system temperature, ?i p = 2 
is the number of polarisations, B = 128 MHz is the observing 
bandwidth, P — 531 ms is the pulsar period and W e q — 10 ms. 
We averaged a total of 202 radio-on and 6 radio-off observations, 
which correspond to total integration times of T — 1174.2 min 
and 33.6 min for the radio-on and -off phases respectively. The 
resultant time-averaged profiles are shown in Fig. [2] We place a 
limit on the mean flux density in a radio-off phase S s ^ 0.022 ± 
0.004 mjy (SNR ~ 3), which is approximately 100 times fainter 
than that of the radio-on phase 5 on = 2.2 ± 0.4 mjy (SNR ~ 
190 0fl Converting these para meters into pseudo-luminosities, us- 
ing ( lLorimer & Kramei]|2005l) 

ii4oo = S1400 d 2 > (2) 

we find 1/1400, off < 2.9 fijy kpc 2 and L1 400, on ~ 0-29 mjy kpc 2 
(where the pulsar distance d ~ 0.36 kpc; lGwinn et alj|l986t) . We 
note that the pseudo-luminosity of PSR B0823+26 in the radio- 
off phase is at least six times fainter than the weak est known radi o 
pulsar PSR J2 144- 3933 (L1400 ~ 20 itJy kpc 2 ; lLorimeilll994h . 
Although this implies that the radio-off phases are consistent with 



4 We note that the val ue quoted here for the radio-on flux density is lower 
than that obtained by lLorimer et a?] <1995t) at 1400 MHz, i.e. 51400 = 
10 ± 2 mjy. We attribute this discrepancy to our longer data set, which 
provides a more rob ust estimate due to the strong scintillation behaviour of 
the source (see, e.g. JWang et alj|2005h . 




Pulse Longitude (deg) 



Figure 2. Average profiles of PSR B0823+26 for the radio-on (top) and 
radio-off (bottom) observations respectively, which are offset for clarity. 
The maximum pulse intensities are normalised to one. The inset plots show 
zoom-ins of the inter-pulse (left) and post-cursor (right) emission in the 
radio-on profile, using the same axes units as the main plot. Note that the 
y-axes of the inset plots represent the intensity offset from one. 

emission cessation, we cannot rule o ut the possibility that they may 
exhibit extremely faint emission (c.f. lEsamdin et al.ll2005t) . 

In Fig. [3] we show an 'activity plot', formed from the 153-d 
one-bit time-series data, which indicates when the pulsar was ob- 
served and whether it was radio-on or -off. The radio emission in 
PSR B0823+26 is clearly observed to undergo modulation over 
variable timescales, which are mu ch shorter than thos e seen for 
PSR B1931+24 (i.e. days to weeks; iKramer et al.l2006h . 

The short-term emission modulation of PSR B0823+26 was 
probed using the data intervals of continuous observations, which 
are detailed in Table [2] These three data sets, shown in Fig. [4] dis- 
play some evidence for 'emission flickering'; that is, rapid changes 
between the radio-on and -off states of emission before a constant 
emission mode is assumed. We note, however, that further high- 
cadence observations are required to confirm this, and the typical 
characteristics of the 'flicker' pulses. 

We also determined the emission phase durations for the ob- 
ject, which are defined as the difference between the start and end 
points (i.e. transition times) of emission phases. To reduce system- 
atic error, the transition times for emission phases are assumed to 
be the mid-points between consecutive radio-on and -off observa- 
tions, that is when the pulsar changes from radio-on to -off and vice 
versa. Fig. [5] shows the result of this analysis. These data clearly 
show that the distribution of radio-on phase durations exhibits a 
broader spread of values compared with that of the radio-off phase. 
Overall, however, there is a bias towards shorter emission phase 
durations, which is in stark contrast to PSR B 193 1+24. The av- 
erage time that PSR B0823+26 exhibited detectable radio emis- 
sion was 1.4 ± 0.4 d. Whereas, the average radio-off timescale was 
0.26 ± 0.04 d. The activity duty cycle (ADC), the percentage of 
time in the radio-on phase, was calculated from the ratio of the to- 
tal radio-on duration to the total observation time. The uncertainty 
in this value was determined from the ratio of the standard error 
in the radio-off time to the mean radio-off time. Subsequently, the 
pulsar is found to be radio-on for 80 ± 10 % of the time. 

We were also interested in determining whether the pulsar ex- 
hibits any systematic t rend in pulse intensity before (or after) a null 
(c.f. P SR B0809+74; [Lvne & Ashworthll 19831 ; Ivan Leeuwen et al.l 
|2002|) . Ideally, this analysis should be performed on single-pulse 
data. As these data were not available, however, we were limited to 
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Figure 3. The sequence of observations of PSR B0823+26 carried out over the 153-d period, denoted by the black lines. The data are separated into three 
continuous N ~ 54 d panels. The times of observation and the times when PSR B0823+26 was radio-on (full-amplitude) and -off (half-amplitude) are shown 
by the extent of the black lines. The times of more intensive observing sessions are shown at M JD ~ 54838, 54855 and 54858. 



54857.7 54857.8 54857.9 54858.0 54858.1 54858.2 54858.3 

Modified Julian Date (day) 

Figure 4. Three sequences of continuous observations of PSR B0823+26 during January 2009. The times when the pulsar was radio-on and -off are shown 
by the extent of the black lines. 
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Figure 5. Histograms showing the durations of time when the pulsar is ob- 
served to be in a specific emission phase, on (top) and off (bottom), with 
inset plots of the radio emission activity for timescales up to 15 h. 



using the 6-min integrated pulse profiles to discern any systematic 
brightness variations in the emission phases. We computed the peak 
flux density for each profile, using the first term in Eqn.Q] assuming 
a constant observing system (i.e. with stable T sys , {3 and G param- 
eters). We find no evidence to suggest any correlation between the 
pulse intensity and pulse integration number preceding (or follow- 
ing) a radio-off phase in these data. In addition, no correlation was 
found between emission phase length and pulse intensity. 



3.2 Periodicity analysis 

In order to elucidate the behaviour of PSR B0823+26, we have per - 
formed weighted wavelet Z-statistic (WWZ) analysis dFosten 19961) 
on the one-bit time-series data. The WWZ algorithm adop ts a mod- 
ified approach to traditional wavelet analysis (see, e.g.. lAddisonl 
120021 for a review of wavelet transforms) in order to counter the 
undesired eff ects of uneven time sampling (e.g. spectral leakage; 
IScargle|[i982l) . It employs a wavelet function which includes a pe- 
riodic, sinusoidal test function, of the form e l "( t_T ), by project- 
ing the data onto a set of sine and cosine trial functions (i.e. the 
waveform). It also utilises a Gaussian wi ndow functio n (weighting 
function of the data) which is defined as (Foster 1996) 



w(u), t) — e 



(3) 



and is centred at time r, with a width defined by the frequency uj 
and tuning constant c. 

As such, the WWZ uses the sinusoidal wavelet to fit the data 
and the sliding window function to weight the data points which, 
in turn, mitigates spectral leakage. We note here that data points to- 
wards the centre of the window in the fit are weighted the heaviest, 



and those near the edges of the window the least. The spectral con- 
tent of a signal is, subsequently, ob tained at times corresponding to 
the centre of the wavelet win dows dBedding et al.l 19981 ; iTempletonl 
l2004l;lTempleton et alj|2005t) . 

For the purpose of our analysis, we chose c = 0.001 so that we 
could strike a balance between frequency and time resolution. The 
resultant WWZ transform of the one-bit time-series data, showing 
the spectral power at successive epochs (or time lags), is displayed 
in Fig. [6] It is clear that the WWZ transform exhibits very spo- 
radic structure. This is thought to result from the effect of irregular 
data sampling on the projection (i.e. local matching) of the WWZ 
wavelet function. During the first ~ 40 d of the data-set, the ob- 
servation sampling is at its greatest, with no gaps and ~ 5 — 10 
observations per day (or more). The three continuous observing in- 
tervals are also included in this date range. Consequently, the data 
obtained over the first ~ 40 d accounts for ~ 62 % of the total. We 
find that the data sampling in this observing period is sufficient for 
the matched wavelet function to return power at several fluctuation 
frequencies, which can be attributed to the locations of the con- 
tinuous observing sessions. This is highlighted by the split in the 
WWZ data, centred around MJD ~ 54948, which shows that the 
locations of dominant spectral features strongly correlate with the 
midpoints of the continuous observing sessions. The most promi- 
nent spectral feature can be seen at approximately 0.44 — 0.36 d _1 , 
corresponding to a period of about 2 — 3 d. As the data sampling 
worsens with time, becoming more irregular (roughly a few obser- 
vations per day) with occasional gaps (one is more than 10 d), there 
becomes a point at which the prominent fluctuation frequencies of 
the system are no longer resolved (i.e. around MJD ~ 54972). 

In order to clarify the significance of these variations, we sim- 
ulated several data-sets using the observed data sampling. We cre- 
ated model data-sets with single periodicities, ranging from 0.25 d 
to 10 d, and analysed them with the WWZ. We also performed 
WWZ analysis on data with randomised activity values. For the 
random data, the spectral power is distributed across the entire 
WWZ plane and no dominant periodicities are recovered. For the 
data with an intrinsic periodicity, however, we find that the funda- 
mental frequency is resolved at virtually all epochal However, the 
fundamental spectral component is also accompanied by a broad 
distribution of power and is substantially less significant after the 
first 40 C0. For simulated periods of the order of a day, we note 
that the resulting distribution of power is similar to the broad 
distribution observed in the real data. This indicates that our re- 
sults are modulated by data- windowing effects. We note that the 
higher period components in the observed data (~ 10 — 30 d) 
are most probably not intrinsic to the pulsar. This is because the 
analysing wavelet has a full window-width (~ 10 cycles) that is 
of the order of the length of the data-set and, therefore, is not suit- 
able to provide m eaningful results about any long-term periodicity 
dTempletonll2004b . These results suggest that there are a number 
of fluctuation frequencies present in the data, which are only re- 
solved during epochs when there is the most frequent observation 
sampling. However, we are cautious to attribute these periodicities 
to the intrinsic variability of the source, due to the influence of the 
observation sampling on the results. 

In order to elucidate the periodicities in the pulsar radio emis- 



5 Around MJD « 54913 — 54922 there is a gap in the observations, 
which often results in a prominent decrease in the WWZ power. 



sents its significance 



I to a time-frequency c 
<Beddingetal.ll998t) . 
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Figure 6. WWZ transform of the 153-d PSR B0823+26 radio emission activity data-set (left) and the corresponding integrated power spectrum (right). During 
the initial ~ 40 d of the data-set, the observation sampling is sufficient to resolve several prominent features. These spectral components are also accompanied 
by a broad distribution of power towards higher frequencies, which are likely analogous to Fourier harmonics. At later times (MJD > 54872), the reduced 
data sampling does not allow insight into the dominant fluctuation frequencies. 



sion, we now consider the error estimation of the WWZ data 
(Fig. |6). We note, however, that such e rror estimatio n is a non- 
trivial procedure. This is highlighted by iFostej dl996h . who state 
that analytic description of errors in a WWZ is very intricate, due 
to the nature of the weighted parametric projection. In addition, the 
assumptions made whilst devising the WWZ statistic practically in- 
validate the formal errors. For example, we assume the null hypoth- 
esis that the data is purely sinusoidal with constant frequency and 
amplitude, plus random noise, which we know is false. Therefore, 
any analytic errors which are calculated are consequently subject 
to these assumptions. As a result, we have applied a more heuristic 
approach to the estimation of WWZ errors. Here we employ two 
different methods, both of which are complementary. 



3.2.1 Confusion limit estimation method 

Following iTempleton et al.l d2005l) . we used the confusion limit 
method to estimate the peak fluctuation frequencies and their maxi- 
mum 1-er uncertainties. Here, the peak fluctuation frequency at a 
given epoch is one which is associated with the greatest WWZ 
power. The maximum l-cr uncertainty associated with this fre- 
quency is approximated by measuring the confusion limit of the 
WWZ spectrum, that is the half-width at half-maximum of the Z- 
statistic, Z{ui, r). The result of this analysis is shown in Fig. [7] We 
find that the average error in the fluctuation frequency from this 
analysis is ~ 11 %. 

The variation in the dominant periodicity, within the first 40 d, 
is above the error limit. However, as the data is clearly modulated 
by data-windowing effects (see above), we believe that these vari- 
ations in the dominant periodicity are most likely governed by the 
observation sampling, rather than any modulation in the pulsar ra- 
dio emission. Through comparing the observed data with the ran- 




1 54840 54860 54880 54900 54920 51910 5 1960 54980 ' 
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Figure 7. The peak fluctuation frequencies (crosses) and periods (open cir- 
cles) from the weighted wavelet Z-transform data. Error bars are l-cr values 



computed using the confusion limit estimation method (fractional errors for 
individual epochs). The average error in the peak fluctuation frequency (pe- 
riod) is ~ 11 %. For the first 40 d, the peak fluctuation period is typically 
~ 2 — 3 d apart from at a couple of epochs (~ 0.8 d), which may repre- 
sent a variation in the intrinsic periodicity. After MJD > 54872, the data 
sampling becomes poorer, and more irregular, which favours the smaller 
(non-physical) peak fluctuation frequencies. 

dom simulated data, we conclude that the source does exhibit some 
(quasi-)periodicity, but which cannot be accurately resolved here 
due to the aforementioned data-windowing effects. 



3.2.2 Data-windowing method 

Uncertainties in the WWZ transform were also estimated by com- 
puting the standard deviation of windowed data. Here, we separated 
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the WWZ data into segments of length T = 14 d, and calculated 
the standard deviation in the peak frequencies (<r„( T )) and peri- 
ods (o>( T )) for each segment accordingly. This provided a mea- 
sure of the modulation in the peak frequency and period over time. 
We note, however, that we only included peak fluctuation periods 
which were below 10 d in this analysis, as these periods are the 
longest which are well represented by the WWZ (see above). 

Table [3] shows the results of this analysis. The median peak 
fluctuation frequencies (periods) for the first three segments of the 
transform data are ~ 0.42 - 0.35 d" 1 (~ 2 - 3 d). We estimated 
the significance of perio dicities within the WWZ data using a boot- 
strapping approach (e.g.lShao. J. & Tu D J 19 95: Politi s. D.N, et al.l 



119991 : IZoubir. A.M. & Iskander. D.rJ|2004|) . Here, we resampled 
the transform data 100 times and, for each resample, determined the 
standard deviation. This allowed us to determine an accurate esti- 
mate for the average standard deviation, or background noise, of the 
transform data from the population of possible values. We assume a 
5 crwwz ~ 80 level as a confident signal detection. Consequently, 
we find that the median-peak WWZ values of the windowed data 
are significant for MJD < 54872. Whereas, those afterwards do 
not meet the cut-off criteria. This further indicates the importance 
of high-time resolution on the WWZ analysis; that is, periodicities 
in the data can only be accurately constrained when sufficient time 
resolution is available. 

3.2.3 WWZ analysis summary 

The above results indicate that PSR B0823+26 may exhibit a num- 
ber of (quasi- (periodic features in its radio emission; the most 
prominent of these being around 2 — 3 d. However, evidence for 
their existence is only obtained during a small portion of the data- 
set (i.e. in the first 40 d), when the observation sampling is at its 
highest. From the simulated data-sets, it also appears that these 
'significant' features may be influenced by spectral leakage (i.e. 
data-windowing effects on the WWZ). In light of this, we postulate 
that PSR B0823+26 does exhibit a fundamental (quasi)-periodicity 
in its radio emission, but stress that further higher-cadence obser- 
vations are required to better characterise this behaviour. 



4 TIMING BEHAVIOUR 



In PSR B 1 93 1 +24, we see systematic variations in its timing resid- 
uals, which manifest as quasi-periodic cu bic structure due to spin- 
down rate variation dKramer et al. I l2006h . To determine whether 
PSR B0823+26 exhibits similar behaviour, we analysed the timing 
measurements obtained from the 153-d AFB data-set. We obtained 
the topocentric TOAs using PSRPROF0 and analysed them using 
the TEMP02 packagtjfl The resulting residuals (observed — pre- 
dicted TOAs) from a best-fit timing model (Table |4]l are shown in 
Fig.® 

We find that there is no significant evidence for any peri- 
odicity in the 153-d AFB timing residuals. However, this is not 
that surprising considering the timescales of emission variation 
(~ hours) will result in an observational bias against such detec- 
tion; the ability to resolve a discrepancy between observed and pre- 
dicted TOAs will be strongly dependent on the length of time a 



7 http://www.jb.man.ac.uk/ pulsar/observing/progs/psrprof.html 

8 A d etailed overview of this timing package is provided by [Hobbs et al J 

Further details and documentation can also be found at 
http://www.atnf.csiro.au/research/pulsar/tempo2/ 



Table 4. The properties of PSR B0823+26 obtained from timing measure- 
ments of the 153-d AFB data-set. Note that the right ascension, declination 
and dispersion measure of the source are held fixed in the fit to these data. 
The standard 1-cr errors are provided in the parentheses after the values, in 
units of the least significant digit. 



Parameter 


Value 


Right Ascension (J2000) 


0H h 26 m 51?489 


Declination (J2000) 


+26° 37' 23 '.'706 


Epoch of frequency (modified Julian day) 


54909.0 


Rotational frequency v (Hz) 


1.884439516337(1) 


Rotational frequency derivative i> (s -2 ) 


-5.997(1) xlO" 15 


Dispersion Measure DM (cm -3 pc) 


19.464 



CO 
"D 



CD 



w) o 




54840 54880 54920 

Modified Julian Date (day) 



54960 



Figure 8. Post-fit timing residuals for PSR B0823+26 from the 153-day 
AFB data-set, after fitting for v and v. There is no apparent cubic structure 
in the timing residuals. 



pulsar assumes a given spin-down rate. Nevertheless, it is possible 
that the object alternates between spin-down rates, consistent with 
the magnetospheric-state changing scenario dLvne et al]|2010l) , but 
we are not sensitive to the variations in the timing residuals. This 
idea will be explored in greater detail in the following sections. 



4.1 Overview of the timing model 

To determine whether the timing measurements of PSR B0823+26 
could be consistent with a variable spin-down model, we developed 
a simulation tool to reproduce the timing behaviour of intermittent 
pulsars. This tool is based on a Monte-Carlo method, whereby we 
define a parameter space of trial radio-on and -off spin-down rates 
(i> on and !>off) to produce simulated timing residuals that are com- 
pared with the observed. The tool also uses an initial rotational fre- 
quency, vq, and average fitted rotational frequency derivative, !> av , 
to simulate the timing behaviour of the object. These additional 
parameters are determined from fitting the pulsar's rotational and 
orbital parameters from the observed data using TEMP02 (see Ta- 
bleH. 

The emission activity of the source is represented by two ana- 
lytic, exponential functions that separately model the observed dis- 
tributions of emission phase durations (i.e. radio-on and -off). The 
pulsar is modelled to sequentially switch between the radio-on and 
-off phases. Therefore, the simulated emission phase durations are 
obtained from the analytic functions by alternately and randomly 
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Table 3. Summary of the results from the data-windowing error analysis. The time span of the WWZ data analysed is denoted by 'MJD range'. The median- 
peak WWZ values, of these segements, and their 1-cr uncertainties are denoted by WWZ max and A (WWZ max) respectively. The peak fluctuation frequen- 
cies, periods and their corresponding 1-cr uncertainities are given by ^fl uc , A(^ uc ), Pfl uc and A(Pq uc ) respectively. 



MJD range 


WWZ max 


A(WWZ max ) 




A(^ fluc )(d- 1 ) 


Pfluc(d) 


A(P fluc )(d) 


54832-54846 


94 


25 


0.35 


0.39 


2.9 


1.0 


54846-54860 


125 


9 


0.42 


0.05 


2.4 


0.3 


54860-54874 


1 13 


32 


0.40 


0.12 


2.5 


4.0 


54874-54888 


63.9 


0.2 


0.0728 


0.0001 


13.79 


0.02 


54888-54902 


62 


1 


0.0728 


0.0001 


13.75 


0.02 


54902-54916 


57 


1 


0.04 


0.02 


28.0 


7.4 


54916-54930 


57.4 


0.3 


0.0357 


0.0000 


28.0 


0.0 


54930-54944 


58.0 


0.1 


0.0357 


0.0000 


28.0 


0.0 


54944-54958 


58.1 


0.1 


0.03570 


0.00004 


28.01 


0.03 


54958-54972 


57.5 


0.3 


0.0356 


0.0000 


28.09 


0.00 


54972-54986 


56.3 


0.5 


0.0356 


0.0000 


28.09 


0.00 



sampling from the possible ranges of values. This process is con- 
tinued until the emission activity of the total observed data dura- 
tion (~ 153 d) is fully covered for each simulation trial, thus pro- 
ducing pseudo-random generated number distributions for both the 
radio-on and -off emission phase durations. This method, in turn, 
provides us with a better general description of the pulsar emis- 
sion activity, due to the finite number of data points we can sample 
from the parent distributions. With the above in mind, the timing 
residuals of each simulation trial represent one possible observed 
outcome, for a given pair of spin-down rates, considering a pulsar 
which exhibits variable emission activity with given distributions 
of switch durations. 

As the spin-down rate of the model pulsar alternates between 
consecutive emission phases, the rotational frequency is updated at 
each integer step in pulse number: 

TOA = - + tret ■ (4) 
V 

The reference time t re f corresponds to the total time elapsed at the 
last step and n is the number of pulses in each emission mode. The 
rotational frequency v is updated using 

V = V Ic f + (i'phase X At) , (5) 

where i/ re f is the reference frequency at the last step, !> pna sc is the 
spin-derivative in the corresponding emission phase and At is the 
duration of the emission mode. 

The simulated TOAs, obtained from these data, are those that 
would be measured at the Solar System Barycentre (Barycentre ref- 
erence frame) and are inclusive of additive white Gaussian noise, 
to simulate instrumental noise. For each simulated TOA, the noise 
signature is calculated from a randomly selected error bar from the 
observed data; a random number is generated from a Gaussian dis- 
tribution with a full-width at half-maximum that is twice the size 
of the error bar. To track pulsar rotation in time, TOAs for radio-on 
and -off phases are calculated. The radio-off phase TOAs are the 
theoretical TOAs which an observer would measure if the pulsar 
was detectable. 

We note here that, as v is evolved over time, the initial value 
for this quantity is equal to the rotational frequency at the start point 
of each simulation trial i/o, which is obtained from the average v 
of the observed data. Therefore, for each combination of i> on and 
!> gf, the simulated value for i/q will be different compared to the 
observed due to the relative contributions of v OVL i0 g. This, in turn, 
results in a systematic offset in the resulting timing residuals. To 



correct for this effect, we fit an average v to each set of simulated 
timing residuals using TEMP02. It is important to note that we 
only consider the radio-on phase TOAs in these fits, so as to ac- 
curately simulate the observations. Following this procedure, we 
then use TEMP02 to provide post-^-fit root-mean-square (here- 
after, simply referred to as RMS) values for each v ori >Q g trial com- 
bination. The optimal combination of v OTL and u a s is then taken as 
that which obtains an average RMS from these trial values which 
is closest to the observed. 



4.2 Simulation analysis 

To calibrate the simulation tool, and test its functionality, it 
was applied to the p rototype intermittent pulsar PSR B1931+24. 
I Kramer et alj J2006j) show that there is clear quasi-periodic, cu- 
bic structure in the timing residuals of this pulsar and fit the data 
to obtain v on = -16.3 ± 0.4 x 10" 15 s~ 2 and i> oS = 
— 10.8 ± 0.2 x 10~ 15 s~ 2 . As a consistency check, we analysed 
the same data-set using a number of trial spin-down rates within 
the stated uncertainties. TOAs obtained between 4 May 2003 
and 9 October 2003 with the Lovell telescope were analysed us- 
ing TEMP02 to determine i/q and z> av , and to create an accu- 
rate ephemeris, to perform the analysis. We find that the simu- 
lation tool reproduces the timing behaviour of PSR B 193 1+24 
well, as shown in Fig. [9] for v on — —16.1 x 10~ 15 s~ 2 and 
z> ff = —10.8 x 10~ 15 s~ 2 . We note here that there are a cou- 
ple of small offsets between the observed and simulated residuals 
(at MJD~ 52820 and ~ 52918), but emphasise that these do not 
reflect on the validity of the simulation tool. We stress that the ex- 
ample shown (Fig. [9]l does not represent a perfectly optimised re- 
sult. It is also important to note that the emission phase durations of 
PSR B 193 1+24 are only known to an accuracy of approximately 
±1 d. This, in turn, will naturally result in a discrepancy between 
the observed and simulated data, if the actual timescales of emis- 
sion are different to those determined by the observation sampling. 

Following the above approach, we obtained vq and i> av for 
PSR B0823+26. We determined the constraints on the combina- 
tions of !> on and z> ff from the pulsar's ADC and average fitted v, 
which by definition is 

Z> a v = (ton X !>on) + (toff X ^off) , (6) 

where t on and t D g are the fractional times that the pulsar is 
in the radio-on and -off emission phases respectively. Accord- 
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52760 52800 52840 52880 52920 

Modified Julian Date (day) 

Figure 9. Observed (crosses) and simulated (open circles) timing resid uals 
of PSR B 193 1+24, using the same data as that in lKramer etal] <2006h but 
with slightly different spin parameters. The results of the simulation are 
consistent with the timing behaviour of the source. 




i> on (irr L v 2 ) 



Figure 10. Residual contour map obtained from simulations of the tim- 
ing behaviour of PSR B0823+26, for given combinations of rotational fre- 
quency derivatives. The plot shows the distribution of RMS values for the 
specified range of on Q ff , with the contours denoting RMS levels for the 
given parameters. The diagonal region in the middle of the plot shows the 
area of acceptable u combinations. The observed timing RMS is 80 ±10 ^ts. 



ingly, the parameter space of spin-derivatives to test was de- 
fined as -6.06 ^ z> on (KT 15 s~ 2 ) sC - 5.99 and 
—6.0 ^ Vofi (10~ 15 s~ 2 ) ^ —5.6 (assuming a small variation in 
the spin-down rate between emission phases). The spin-derivatives 
are inherently constrained by !> av ~ — 5.997 x 10 -15 s -2 , so 
that neither i> on or v g can cross this boundary. For each combina- 
tion of input parameters, the RMS was calculated for 2000 itera- 
tions of the simulation tool using a resolution of 2 x 10~ 18 s~ 2 
in Von, off- These data were then averaged and compared with the 
RMS and measured uncertainty from the observed data. The un- 
certainty was calculated using TEMP02 to fit the observed data 
across individual segments. Here, the average and standard devia- 
tion of the sample population of observed RMS values were com- 
puted, and then used to obtain the fractional error and measured un- 
certainty in the RMS. Assuming a 3-cr error in the observed RMS, 
we obtain a constraint of RMS = 80 ± 10 /is on the simulated 
results. 

The result of this analysis can be seen in Fig.|T0j The parame- 
ter space of i/ on — v g where the RMS converges on the observed 
value is seen to be distributed, within the observational errors, to- 
wards the central diagonal portion of the plot. This can be attributed 
to Eqn. [6] which shows that i> av of the simulated data, and hence 
the RMS by logical progression, will converge with the observed 
for z> on oc l/z> ft- Consequently, we find that the maximum in- 
crease from !> D ff to !> on (A;> max ) which satisfies the model criteria 
is approximately 6 %. 

We also note that as the emission phase durations are sampled 
randomly from the model distributions, there exists a variance in 
the computed average RMS for each combination of input param- 
eters. As such, we find that the results of the simulation tool are 
less variable for v values which are closer to the observed average; 
the variance in the results from the simulated residuals has to be 
smaller for the simulation to be consistent with the observations. 
Despite the fact that the simulation tool can be applied to any inter- 
mittent pulsar, it is ideally suited to modelling those with small Az> 
as least-squares fitting methods, such as the method implemented 
bv lKramer etai] feOOfjl) . will be more efficient for larger Ai>. That 
is, the simulation tool can provide an upper limit to Ai>, which is 
independent of curve-fitting analyses which are unfeasible when 
analysing timing signatures with negligible cubic structure. 



4.3 Switching timescale dependency 

Further to the analysis of the simulated data, the effects of alter- 
ing the emission phase durations and ADC were investigated. We 
were interested in determining whether our results were subject to 
systematic effects and how they were affected by modelling the 
emission behaviour differently. 

Firstly, we investigated the effect of over or under-estimating 
the observed emission phase durations. We added systematic er- 
rors to the observed durations, t b S! with a magnitude from 
% to 100 % of t bs in steps of 5 %, with the polarity of the er- 
ror (for each emission mode) being determined by a random num- 
ber generator. As the percentage error was increased, the average 
RMS over 2000 iterations of the simulation tool was found to in- 
crease also. That is, as the variance in the emission phase dura- 
tions increases, so do the fluctuations in the timing residuals due to 
the model pulsar assuming a spin-down rate for longer or shorter 
than is observed. We find that the tolerance of the RMS to these 
systematic errors increases towards smaller variations in the spin- 
down rate. This is particularly true for the combinations of spin- 
parameters which were comparable with the average spin-down 
rate; the results of the simulation remained consistent with the ob- 
served even for an error of ± 100 % in the emission phase dura- 
tions. 

To confirm the validity of the maximum allowed change in 
spin-down rate, Ai> max ~ 6 %, we also performed a simulation 
trial using the observed switching times without any systematic 
errors. Here, we used the corresponding values for v on i0 g which 
resulted in the maximum spin-down rate variation. From this anal- 
ysis, we obtain an RMS which is greater than that observed. How- 
ever, this result is only derived from one trial and, as such, is not 
expected to be as consistent as averaging over many simulated re- 
sults. In addition, the value quoted for A!> max inherently serves as 
an upper limit. Therefore, it is possible that the pulsar may undergo 
more modest variations in its spin-down rate; simulations incorpo- 
rating these spin-parameters produce results that are more consis- 
tent with observations. We note, however, that systematic errors in 
the observed emission durations will have a significant effect on 
the Ai> max limit obtained. As a result, we argue that timing studies 
of intermittent pulsars will be more conclusive for better sampled 
data-sets. 
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We also tested our assumption that the emission activity of 
PSR B0823+26 can be modelled by randomly sampling two distri- 
butions of emission phase durations. Here, we compare the results 
of separate trials which used the analytic and observed (both or- 
dered and randomised) sequences of emission phase durations. We 
find that the average RMS values from these trials are all consistent. 
Therefore, we conclude that the assumptions of our simulation tool 
are valid. 

In addition, we examined the dependency of the results on the 
ADC assumed. In this analysis, the sequences of radio-on and -off 
emission phase durations were combined together and randomised, 
so that no distinction was made between the different distributions. 
This, in turn, resulted in an ADC which no longer coincided with 
the observed value, that is ~ 50 %. A number of v combinations 
were used in this analysis, from A;> m i n to A;> max . However, no 
RMS values within the observed 80 ± 10 /is could be obtained (ex- 
cept for !> on ~ Voft )- This indicates that the form of the distribution 
which is used to model the emission phase durations is fundamental 
to the simulation procedure. 

To investigate this further, and to determine how unique the 
A^ max solution is, the analytic distributions of emission phase du- 
rations were also directly altered. Here, simulations were carried 
out where the parameters of the distributions were modified to ob- 
tain two test cases, for ~ 5 % decrease and increase in the ADC. 
These simulations returned u on . Q ff matches that corresponded best 
with the ADC. That is, for an approximately 5 % decrease in the 
ADC, the maximum change in the spin-down rate reduces, and vice 
versa. For a lower percentage radio-on time, Ai> is required to be 
smaller to obtain !> av , and the opposite for a higher percentage. As 
a result, we obtain limits on the maximum change in spin-down 
rate 4.5 % < A z> max < 7.6 %, if we assume a maximum error of 
~ 5 % in the modelling of the random number distributions. Conse- 
quently, we see that alteration to the ADC has a significant effect on 
the results; if the initial parameters used are wrong, the simulation 
tool will not converge on an optimum combination of spin-down 
rates. Despite the incorporation of a 5% uncertainty in the ADC, 
it is apparent that the upper limits on the spin-down rate variation 
are still small com pared with that of PSR B 193 1+24 (~ 50%; 
iKramer et ai]|2006h . 



5 DISCUSSION 

5.1 Nulling timescales in PSR B0823+26 

While our data represent the most comprehensive study of the long- 
term nulling behaviour of PSR B0823+26, we suspect that the 
majority of the inferred nulling timescales are influenced by data- 
windowing effects. That is, we predict that the emission phase dura- 
tions are often overestimated due to insufficient time coverage. This 
is strongly supported by the activity plots (Fig.[3]— |4j which, when 
compared, clearly show that the pulsar exhibits a much shorter 
modulation timescale (~ minutes to hours) than can be resolved 
in the typical observing cadence (~ days). This effect also influ- 
ences the result of the WWZ analysis (Fig. [6}, due to the lack of 
any significant spectral features at times of poor observation ca- 
dence. Therefore, we cannot explicitly state the typical modulation 
timescale(s) intrinsic to PSR B0823+26; it is likely that the sparser 
sampling of the typical observations may have resulted in missing 
numerous radio-on and -off states. 

Although the majority of our data is subject to the data- 
windowing effects mentioned above, the data from the continuous 



observing runs is not. Therefore, we can use these data to place 
confident limits on the short-term variability of PSR B0823+26. 
In light of this, we conclude that the source undergoes nulls over 
timescales of minutes to hours (up to at least 5 h), but do not 
rule out the possibility of longer (~ 1 d) null phases, which re- 
quire further investigation. Remarkably, in these observing runs, 
we also find evidence for highly irregular short-term (< 1 h) mod- 
ulation. That is, the object appears to exhibit periods of sponta- 
neous bursts of emission and nulls, or emission flickering, before 
a constant emission phase is assumed. If this is typical for neutron 
stars, it may be inferred that pulsars which exhibit intermittency, 
e.g. normal nulling puls ars and rotating radio transients (RRATs; 
iMcLaughlin et al.l 20061 ). may undergo a radio emission 'ignition 
phase' before attaining the necessary criteria for stable emissioifl 
However, this behaviour could just be intrinsic to PSR B0823+26 
and, therefore, requires the study of other sources to associate it to 
the pulsar population as a whole. Evidently, this would be comple- 
mented by further investigation into the typical nulling behaviour of 
PSR B0823+26, which may provide insight into this phenomenon 
and the mechanism(s) responsible for emission modulation in pul- 
sars. 

Taking in mind that nulling pulsars exhibit null durations of 
1 — 10 pulse period s up to approximately a year (PSR J 184 1—0500; 
ICamilo et af]|2012l) . PSR B0823+26 may be considered to be a 
bridge between the different 'types' of nulling pulsars. That is, 
the object may be located inbetween conventional nulling pulsars 
and longer-te rm intermittent pulsars on a 'nulling continuum' scale 
jKeandlioiOl) . A question to now ask then, is if these sources all 
show different cessation timescales will all of them exhibit Ai> and, 
if so, how will it manifest? This will be investigated in greater detail 
below. 



5.2 Evidence and implications of Ai> in PSR B0823+26 

For the first time, we have presented a simulation tool which 
can model the timing behaviour of an intermittent pulsar and 
obtain an upper limit on its spin-down rate variation. Through 
simulating the rotational behaviour of PSR B0823+26, we find 
Ai> max ~ 6% between emission phases, which corresponds to 
i> on = -6.05 x 10~ 15 s" 2 and v oS = -5.702 x 10" 15 s -2 . 
However, we cannot discount a scenario where there is no varia- 
tion in spin-down rate between the different radio emission phases. 
Keeping this in mind, we can still try to estimate the implied cur- 
rent flow in the pulsar magnetosphere for b oth the radio-on and 
-off states using the above values. We follow IKramer et all J2006h 
and consider the simplest possible emission model. That is, we as- 
sume that in the radio-off state the pulsar spins down by a mech- 
anism that does not involve substantial particle ejection (e.g., via 
magnetic dipole radiation if the pulsar is in vacuum). Whereas, in 
the radio-on state, we assume that spin-down rate is enhanced by 
a torque from the current of an additional plasma outflow. Modi- 
fications of this simple assumption are possible (e.g. considering 
the spin- down contribut ion of a plasma-filled close-field line re- 
gion, see iLi et alj[2012al) . but for our purposes the simplest model 
is sufficient. Hence, with the assumption that the increase of the 
spin-down rate is purely due to the torque of the charged plasma 
additionally existing in the radio-on state, over the radio-off state 
vacuum spin-down, one derives an estimate for the charge density 



9 In this context, RRATs may never attain the criteria for stable emission. 
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of this wind component: 

31 Ave 2 



0.0024 Cm" 



(7) 



1.26 x 10 21 i? 6 V^ff 
where c is the speed of light in a vacuum, R is the radius of the pul- 
sar (taken to be 10 6 cm) and I is the moment of inertia (assumed 
to be 10 38 kg m 2 ). We note that this equation also assumes that the 
object is an orthogonal rotator, i.e. the source has a magnetic incli- 
nation angle of 90° . While this assumption is almost certainly not 
met for most sources, this seems to be indeed a good approxim a- 
tion for PSR B0823+26 (a ~ 86°; lEverett & Weisberdl200llPl 
Comparing the above estimate with the Goldreich- Julian density, 
i.e. the screening density in a simple dipolar spin-down model (e.g. 
lLorimer & Kramej|2005b 



PGJ = 



0.0203 Cm" 



(8) 



where B s is the surface magnetic field strength of the pulsar, we at 
least obtain a good indication that only a small amount (~ 12 %) 
of the total available charge in the wind component is contribut- 
ing to the overall spin-down. This indicates that the open field line 
region of the pulsar magnetosphere may never become entirely de- 
pleted of charge during a radio-off phase (c.f. Ili et al.|[2012al) or, 
alternatively, that the spin-down of the pulsar is dominated by the 
non-wind contribution (i.e. magnetic dipole radiation). 

For PSR B 193 1 +24 it is suggested that the radio emission ces- 
sation results from re-configuration of the global magnetospheric 
charge distribution, and that the spin-down rate variation is a sig- 
nature of this phenomenon i Kramer et al.l 120061 ; iTimokhinl |2010| : 
iLvne et alj[201Cl : lLi et alj|2012al lbl). Certainly, i n several other pul - 
sars, we see evidence to support this theory jLvneetalJboiOl) . 
However, it is not entirely clear what manner of emission modula- 
tion will be ac companied by a va riable spin-down rate; in the study 
performed by ILvne et aT] | |2010|) . pulsars with similar changes in 
spin-down rate were observed to undergo different magnitudes of 
variation in their pulse shape. It is interesting, therefore, to find that 
PSR B0823+26 exhibits the same breakdown in radio emission 
production (or detectability) as PSR B 193 1+24, without the need 
for a large change in the magnetospheric currents. This implies that 
the mechanism which produces radio emission is highly sensitive 
to even the smallest changes in the magnetosphere and that mode- 
changing and nulling are closely related. From this, we infer that 
there could be a significant number of pulsars that exhibit regulated 
changes in their spin-down rate and emission, which contribute to 
'timing noise', and that we are not yet aware of them due to small 
fractional changes in v or short timescale variations. 



5.3 Detection limits on Ai> 

While a large number of pulsars are thought to exh ibit spin-down 
rate variation faobbs et aUl20ld ; Evne et alJPioTol) only 20 are 
known to exhibit any discernible change in ii (e.g. I Kramer et al.l 
l2006l : lLvne et al.ll2010r . lCamilo et alj|2012h . With this in mind, we 
sought to determine our sensitivity to Av in neutron stars and, 
hence, whether we can expect to detect this effect in all mode- 
changing and nulling pulsars. 

We begin by considering a simple model of a dual-!> nulling 



10 In fact, we propose that deviation from a simple spin-down model could 
in principle be used to infer the magnetic inclination angle of a source in an 
independent fashion. 
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Figure 11. The expected change in rotational frequency, A^ psr , due to an 
average spin-down rate deviation, |Ai> av |, that is assumed for a time T. 
Overlaid are lines denoting the possible detection requirements of a vari- 
able spin-down rate in PSR B0823+26, that is (from left to right) Ai/ mot j, 
5 A^ mod and 10 A^ mod . 



pulsar, which alternately exhibits !> on and v q in its corresponding 
emission phases, for a given length of time Ati. If the data is fitted 
with a st andard timing model with a single value of z> (!> av = v plcd ; 
see, e.g- lManchester & Tavlodl977h . we expect the total difference 
between the observed and predicted spin-down rates, in each emis- 
sion phase, to be Az>j = !> bs — ^prcd- Extending this model to N 
emission phases, we would expect the total discrepancy between 
the observed and predicted v to be 

A^p S r = Au on Aton, 1 + Al> oS AiofF 1 + Ai'on Aion, 2 H 



E 



Ai>i AU 



(9) 



In light of this, it will be possible to detect the presence of a 
variable ii in the timing residuals of a pulsar when the condition 
A^p ar 2> Ai^ mo d is satisfied (where Az/ mo d is the uncertainty 
in the model v over a time span T). With the above in mind, we 
computed the expected Av psl over the course of three subsequent 
emission phases (of total length T), assuming an average deviation 
in spin-down rate |Ai> av | = \{v ohs — !> pre d)|, for a wide range of 
parameters (see Fig. 11 It. We note that, while these data do not pre- 
cisely model Ai/p Sr , they do offer useful estimates that can be used 
to predict the detection limits on Ai> in pulsars. 

Assuming PSR B0823+26 exhibits a maximum modulation 
timescale At ~ 1 d (hence T ~ 3 d), and |A!> av | ~ 9.7 x 



10" 



s , we would expect Ai^ p 



2.5 x 10" 



Hz. This is 



approximately 2000 times smaller than Ai/ mo( j from our timing 
measurements (~ 5 x 10~ 8 Hz), which can clearly explain why 
we do not find direct evidence for multiple ii values. Comparing 
this with PSR B1931+24, we find that Av psr /Av mod ~ 3, for 
|Ai> av | = 8.8 x 10~ 16 s~ 2 and (T) ~ 40 d, does result in sig- 
nificant detection of a variable v. This indicates that our method of 
calculation for Au psl is robust, and that it is not large enough in 
PSR B0823+26 to facilitate the detection of a variable v. 

Moreover, considering that the maximum null duration re- 
solved in PSR B0823+26 (during the continuous observing runs) 
is about 5 h, it is likely that Ah> psl is overestimated by a factor 
of ~ 5. Therefore, it is not unsurprising that A^ psr is not signifi- 
cant enough to require the incorporation of another v in the timing 
model. In light of this result, it is extremely likely that Ai> will not 
be discernible in a significant proportion of nulling pulsars. 
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With the above in mind, it is interesting to determine the typ- 
ical properties of a pulsar which can facilitate the detection of Au. 
If we assume a typical pulsar has !> av = —4.46 x 10~ lj s -2 , 
|Ai> av | ~ 4.5 x 10~ 16 s~ 2 (i.e. ~ 10% of !> av ) and a detec- 
tion threshold A!/ psr /Af mo d ~ 5, then we should expect to detect 
spin-down rate variation in pulsars which exhibit T > 26 d and 
have a timing precision Ai/ mo d ~ 10 -9 Hz or better. The majority 
of mode-changing and nulling puls ars, however , do no t exhibit such 
long modulation timescales (e.g. IWang et al.l 120071: iKeane et al .1 



2O10t iBurke-Spolaor & Bailed |201Ct IBurke-Spolaor et all boil! 
Keane et alj|201 ll) . Accordingly, we argue that it will be extremely 
difficult to determine whether v switching is a global phenomenon 
in these sources, especially those with very short modulation 
timescales (i.e. < hours). Clearly, more observations of similar ob- 
jects are required before a consensus is reached on their common 
properties, and how PSR B0823+26 fits into their framework. 



large changes in v. Ultimately, further (multi-frequency) study of 
PSR B0823+26, and other objects like it, should reveal significant 
information about the nature of normal nulling pulsars, RRATs 
and longer-term intermittent pulsars. This, in turn, should offer the 
possibility to determine how all these objects are related (if they 
are), what their typical characteristics are and what mechanism is 
governing their irregular behaviour. 
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6 CONCLUSIONS 

We have found evidence to suggest that PSR B0823+26 exhibits a 
broad distribution of nulling timescales i.e. < minutes up to several 
hours or more. Although longer duration nulls (~ day) have been 
observed in our data, we cannot rule out the presence of short-term 
variations < minutes in the radio-off phases which we are insensi- 
tive to; radio emission flickering, which was observed in the con- 
tinuous observing runs as 'pre-ignition' radio pulses, could have 
occurred during times of low observation cadence and, hence, been 
missed. As such, long (~ hours) single-pulse observations of this 
source are required to confirm its overall nulling fraction, emission 
intermittency timescales and, therefore, its periodicity; that is if the 
source does not exhibit random emission fluctuations. These data 
should also provide further information as to why nulling (a.k.a. 
intermittent) pulsars might undergo emission flickering before they 
assume a stable magnetospheric state. 

A simulation tool was developed to reproduce the timing be- 
haviour observed in intermittent radio pulsars. This tool provides 
evidence for an upper limit of approximately 6 % on the variation 
in the spin-down rate of PSR B0823+26. However, we cannot rule 
out a scenario where the pulsar retains a constant, single v. Nev- 
ertheless, the low limit on Ai> max , compared with that observed in 
PSR B 193 1+24 (~ 50 %; lKramer et alj2006h . suggests the impor- 
tance of small changes to the global charge distribution of a pulsar 
magnetosphere, which might easily perturb the production or de- 
tectability of radio emission. With this in mind, PSR B0823+26 
could provide a link between conventional nulling and more ex- 
treme pulsar intermittency, due to the long timescale radio emission 
modulation and small change in ii between emission phases. 

Despite the somewhat unique emission geometry of 
PSR B0823+26 which facilitates its detection as an inter-pulse 
pulsar, no morphological clue has been provided to explain its 
intermittent behaviour. This is primarily due to a dearth of studies 
which attempt to correlate pulsar emission geometry with pulse 
intensity modulation, as well as the ordinary location of the object 
in the P — P diagram which, ultimately, make it extremely difficult 
to discern the general properties (e.g. P, P and a) of the source 
from the bulk of the normal pulsar population. It is interesting to 
note, however, that the object undergoes intermittent behaviour 
without the need for large changes i n spin-down rate, contrary 
to that predicted by iLi et al.l j2012al) (z> n/j> ff ~ 120 % for 
a ~ 90°). In this context, it may only be pulsars which exhibit 
long modulation timescales (> days to weeks) that experience such 
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